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Abstract. In this paper, we consider the S/RC (Susceptible-Infected-Recovered-Cross immune) epidemic 
model. First the non-negative solution of the SIRC model in fractional order is presented. Then the multi-step 
generalized differential transform method (MSGDTM) is employed to compute an approximation to the 
solution of the model of fractional order. The obtained results are compared with the results by forth order 
Runge-Kutta method and nonstandard numerical method in the integer form. Finally, vve present some 
numerical results. 


1. Introduction 


Infectious agents have had decisive influence on the history of mankind. In fourteenth 
century Black Death has taken lives of about one third of Europes population [1]. The first 
major epidemic in the USA was Yellow Fever epidemic in Philadelphia in 1793, in which 
about 5,000 human died out of total 50,000 population. Epidemic models are used to 
understand the spread of infectious diseases in population 11, 2]. The practical use of 
epidemic models must rely heavily on the realism put into the models. This does not mean 
that a reasonable model can include all possible effects but rather incorporate the 
mechanisms in the simplest possible fashion so as to maintain major components that 
influence disease propagation. Great care should be taken before epidemic models are used 
for prediction of real phenomena [3]. However, even simple models should, and often pose 
important questions about the underlying mechanisms of infection spread and possible 
means of control of the disease or epidemic. Kermack and McKendrick first time introduced 
an epidemic model see for more detail [2]. These papers have had a major influence on the 
development of mathematical models to capture the spread of different diseases and are 
still useful in many epidemic situations. The models presented in these papers laid out a 
foundation for modeling infectious diseases in constant population. It means that there is 
no birth or death from infection. Kermack and McKendrick [2] in their first paper start with 
the assumption that all members of the community are initially equally susceptible to the 
disease, and that a complete immunity is conferred after the infection. The population is 
divided into three distinct classes: S represents healthy individuals who can catch the 
disease; / represents the number of infected individuals; R represents the recovered 
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(removed) individuals who have had the disease and are now immune to the infection (or 
removed from further propagation of the disease by some other means). These 
subdivisions of the population are called compartments. To formulate in regular order, the 
individual goes through consecutive states like $= I— R, 

such models are often called the SIR models. 

Two epidemic models: the susceptible-infected-susceptible (S/S) model and the susceptible 
infected-recovered ( SIR) model are commonly used for studying the spreads of epidemics 
[3,4]. In the SIS epidemic model a recovered individual can be infected again while in the 
STR model which assume recovered individuals have lifelong immunity to the disease and 
this difference makes them suitable for different kinds of infectious diseases. For instance, 
childhood diseases in which individuals can have long-lasting immunity, either naturally or 
from vaccination, are appropriate for S/R model. While for viruses transmitting infection, it 
is more reasonable to use S/S model. Several researchers considered both S/S and S/R 
epidemic models and presented different epidemics in different regions all over the world 
[3-5]. Hethcote [6], presented the interaction of susceptible S(t), infected /(t) and 
recovered R(t) individuals is given by: 











dS (t) 
= N(t)- uS(t) - BS(t)IC(t), 
dt 
dI(t) 
= BS(t)l(t)-(y + u)l(t), 
dt 
dR(t) 
= yI(t)- uR(t). 
dt 


Here M) is the total population, g is the interaction rate of infection, ~ is the death rate 
and ; is the recovery rate. Zaleta and Henandez [7] considered a simple two dimensional 
SIS model with vaccination showing backward bifurcation. Casagrandi et al. [9] presented 
the SIRC (Susceptible-Infected-Recovered-Cross-immune) model. This compartment (C) 
presents an intermediate state between the fully susceptible ($) and the fully protected (R) 
one. For numerical solutions Jodar et al. [9] presented the nonstandard finite difference 
schemes. Also Samanta [10] extended the work of Casagrandi et al. [9] for time dependent 
population size and distributed time delay. But all these work has been done in the integer 
order differential eguations. 

Novvadays several researchers work on the fractional order differential equations because 
of best presentation of many phenomena. Fractional calculus represents a generalization of 
the ordinary differentiation and integration to non-integer and complex order. Fractional 
calculus is used to established new models in many fields not only in mathematics see for 
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example [11-16]. For this purpose Shahed and Alseadi [17] developed a fractional SIRC 
model. In their work they presented a detailed analysis for the asymptotic stability of 
disease-free and positive fixed point. 

In this paper, we consider an SIRC model. First we show the positive solution of SIRC 
model in fractional order. Then we use the multi-step generalized differential transform 
method to approximate the numerical solution. Finally we compare our numerical results 
with nonstandard numerical method and forth order Runge-Kutta method. 

This paper is organized as: In Section 2, we present formulation of the model with some 
basic definitions and notations related to this work. In Section 3, we show the non-negative 
solution and unigueness of the model. In Section 4, the multi-step generalized differential 
transform method (MSGDTM) is applied to the model. In Section 5, the numerical 
simulations are presented graphically. Finally, we give conclusion. 


2. Formulation of Model with Preliminaries 


Here, we consider the model taking by M. El-Shahed etal. [17]. 














dS (t) 
= KO — S(t)) — BI(t)S(t) + ye (7), 
dt 
dl(t) 
= PI(t)S(t) + opC(t)I(t)- (u + 0)I(t), 
dt 
(1) 
dR (t) 
=(1-0o)pPC(r)I(t) + 0It)- (u +0)R(t), 
dt 
dC (t) 
= ORG) —- pPC(t)I(t)-— (u + Y)C 2). 
dt 


With S(0) = Ban Z(0) = Z; R(0) = Re, CO) = Co: 


I , SC l I ai 
Here g is the contact rate of infection, y is the cross-immune period, 9 is the 


, , , -1, , , , , 
infectious period, 6 is the total immune period and s is the fraction of the exposed 
cross-immune individuals who are recruited in a unit time into the infective subpopulation. 
The total population N (t) = S(t) + I(t) + R(t) + C(t), 


so we obtain by adding all equations of the system (1) 
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dN (t) 
sab): (2) 





dt 
Now we introduced fractional order to the system (1) which is consisting of ordinary 


differential eguations. The new system is described by the following set of fractional order 


differential eguations: 


D, S(t) = HO SH) - BSH) + YC (1), (3) 
D, I(t)= PSIG) + GBC(t)I(t) - (u + 0)1(1), (4) 
D, Eis 0G) £C(1)I(t) + OIO) = (y + 6)R(7), (5) 
D, C(t) = 6R(t)- BC(0)I(t)- (+ y)C (t), (6) 
D, NG) u- uN (t): (7) 


Here we used the Caputo sense fractional derivative D’. 


Now we give some basic definitions related to this work and fractional calculus [11-16]. 


Definition A function f(x)(x > 0) is said to be in the space C (a e R) ifit can be written 
as for some p > x Where f(x) Is continuous in10,c), and it is said to be in spacec” if 


pO OC MEN: 


Definition The Riemann-Liouville integral operator ofordera>0 with a >0 is defined as 


X a -1 
fa (x =f) f(t)dt, x> a, (8) 





a 
(1, faye 
I (a) 


0 
Ja JA) = f(x). (9) 


Properties of the above operator can be found in [11]. We only need the following: 


Definition For fe Coi and for a,p > O,a>2 O,ce R and y > —1, we have 


a B p a a +p 
FT OSG A O AU), (10) 


yta 
a y X 
JX = B ,„(4,y +1), (11) 
l(a) x 
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where B (a,y +1) Is incompilete beta function which is defined as 


B_(a,y +1) = ft d-t)" dt, (12) 


Ge = k 
i open) 0 e, 
kso V (a +k +1) 


S 
A 
|! 


(13) 
Definition The Caputo fractional derivative of f(x) oforder a > 0 With a 20 Is defined as 


1 x (m ) 
(m) (t) 
)(x) = ss - dt, (14) 


a +1-m 


(Deg. V 
r(m sea) Kat) 


for m -1<a <m,me N,x2 a, f(x) € C 
The fractional derivative was investigated by many authors, form-1<a<m,f(x)e cC” 


and a > -1, we have 


a a m m GE (k) (za) 
JP J MAH 4 D fla)j= f(x)- 2 f GE (15) 


For mathematical properties of fractional derivatives and integrals one can consult the 
mentioned references. 


3. Non-negative solutions 


Let R?={X eR’ : X >0) and x (t) = (S(t), IC) RG)C (£), NC). For the proof of the 
theorem about non-negative solutions we shall need the following Lemma [16]: 
Lemma (Generalized Mean Value Theorem) Let f(x) Cla,bj and D'f(x)e Cļa,b] for 


O<a<l. 


Then we have, 


1 a a 
D f(Eé)(x - a) (16) 





FARA) = J LE) + 
I (a) 


with OS šā for all xe (a,b]. 


Remark Suppose f(x) € C[a,b] and D'f(x)€ C a,b] for 0<a <1. Itis clear from the 
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above Lemma that if D* f(x)20, for all xe (0,b), then the function f is non-decreasing, 


and if D* f(x)<0, forall xe (0,b), then the function f is non-increasing. 


Theorem There is a unique solution for the initial value problem given by (3)-(7), and the 


solution remains in R”. 


Proof The existence and uniqueness of the solution of (3)-(7), in (0,«) can be obtained 
from [16, Theorem 3.1 and Remark 3.2]. We need to show that the domain R” is positively 


invariant. Since 


D°R|_,=(l-o)BCI + O12 0, 


a 
D C|. = 6R >0, 


On each hyper-plane bounding the nonnegative orthant, the vector field points into R”. 


4, Multi-step generalized differential transform 
method. 


We applying the multi-step generalized differential transform method to fined the 
approximate solution of equations (3)-(7), which gives an accurate solution over a longer 
time frame as compared to the standard generalized differential transform method. Taking 
the differential transform of equations (3)-(7) vvith respect to time vve obtain, 
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(ak +1) k 
e (xa Je Shina ye). 
F((ak +1) + 1) SCH 
T(ak +1) k M 
I(k +1) = qi S(k-i)I(i)+ op X C(k-i)I(i) — (u + DI), 
F((ak +1) + 1) E ge 
(ak +1) k 
R(k +1) = lu —0)B% C(k-i)I(i) + 0I(k) - (u + 6)R(k)), (17) 
F((ak +1) + 1) lā 
T(ak +1) k 
Čiksi)s (8 R (k) dēles PE), 
T((ak +1) + 1) Se? 
(ak +1) 
N (k +1) = ———(4 - uN(k)). 


T((ak +1) +1) 
Here S(k), I(k), R(k),C(k) and _W(k) are the differential transformation of 


S(t), I(t), R(t), C(t) and Nr). The differential transform of the initial conditions are 


S@)=S,, 1(0)=1,,R(0)=R,,C(0)=C, and N(O)=N,. 
In view of the differential inverse transform, the differential transform series solution for 


the system can be obtained as 


E SOG 

| K ak 

| I(t)= 2 I(k)t 

| R(t)= > Rie (18) 
| K ok 

| Ct) = } C(k)t 

[vo Dd", 


Now according to the multi-step generalized differential transform method the series 
solution for the equations (3)-(7) is suggested by 
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| > Sy, te [0,t,] 
| ti ak 
| QI) ei 
Sal t (19) 
| , 
| 
= ak 
E KEE EE WE EE 
| 5 On, te [0,t] 
| K a k 
| e i) ER 
gas | (20) 
| , 
| 
e a k 
E E E 
| Rus, te [0,t] 
EE Ee El 
vw | (21) 
| | 
| 
e a k 
E EE A) EE ās I 


Dei”, eo 


| 
ll, Në 
| Z CODE -= t) LN 

cl | (22) 
| | 

| 

(ze, is ela 


-1 M -1 M 
k=0 
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WE: N GOT, te 10,1,) 
| K ak 
i NO Kë 
"= : (23) 
| | 
| 
e N (kat, EE 


Here S (k), I,(k), R,(k),C (k) and N (k) for j-1,2,..,M Satisfy the following 


recurrence relations 








T(ak+1) k 
S (k+1)= —(ua- S (k)) - PÈS (k-i)I;(i) + re ,(k)), 
T((ak +1) +1) RS 
T(ak+1) k k 
TI (DI (+2 KE (ksi) TG) = (urā oi) 
F((ak +1) +1) E Se 
T(ak+1) k 
R (K-41) = Juegt C j(k=i)1 (i) +61 (k) Ga + SDR (k) ), (24) 
F((ak +1) +1) a 
T(ak+1) k 
C (k +1) = (6R,(k)- PĒC R-DIG) - (a +C), 
T((ak +1) +1) SE 
T(ak +1) 
N ,(k +1) = ——( V 2 uN ,(k)). 


T((ak +1) +1) 


(0), R, DO 


With the initial conditions S,(0)= S, ,(0), I,(0) = I, | 


C ,(0) = C, ,(0), and N (0) = N, ,(0). Finally, we start with initial conditions 


S(0) = S). I(0) = I), R(0) = Ry, C(O) = C, and N(0) = N}, and use the recurrence 
relation given in the above system, we can obtained the multi-step generalized differential 


transform solution given in (19)-(23). 
5. Numerical Method and Simulation 
We solve analytically the equations (3)-(7) with initial conditions by using the multi-step 


generalized differential transform method (MSGDTM). VVe also numerically solve by 
nonstandard numerical method and forth-order Runge-Kutta method for numerical results. 
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For numerical simulation we use a set of parameters given in Table 1. To demonstrate the 
effectiveness of proposed algorithm as an approximate tool for solving the nonlinear 
system of fractional differential equations (3)-(7) for large time ¢ we apply this algorithm 
on the interval [0-30]. 


Table 1: Parameter values for the numerical simulation 


Notation Parameter description Value 
u Natural death rate 0.02 

B Contact rate of infection 0.03 

y Cross-immune period 0.05 

6 Infectious period 0.07 

d Total immune period 0.01 
o Fraction of the exposed cross-immune 


individuals who are recruited in a unit 
time into the infective subpopulation 0.08 


From the graphical results in Figs. 1-5, it can be seen that the results obtained using the 
multi-step generalized differential transform method match the results of the classical 
Runge-Kutta method very well, which implies that the presented method can predict the 
behavior of these variables accurately for the region under consideration. 

Figs. 6-10 show the approximate solutions for S(t), E(t), foi R(t) and N (t) obtained for 


different values of o using the multi-step generalized differential transform method. 
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shows the susceptible individuals 


* MSDF 
—-—--RK4 
— — — ODE45 


susceptible individuals 





0 5 10 15 20 25 30 
Time (day) 


Fig 1. Shows the susceptible individuals. 


Shows the infected individuals 


+ NSDF 
—:—:- RK4 
— — — ODE45 


Infected individuals 





0 5 10 15 20 25 30 
Time (day) 


Fig 2. Shows the infected individuals 
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Shows the recovered individuals 


+ NSDF 
—-—--RK4 
— — — ODE45 


Recovered individuals 





E 5 10 15 ER 25 30 
Time day) 


Fig. 3. Shows the recovered individuals 


Shows the Croass-immune individuals 


+ NSDF 
—-—--RK4 
— — — ODE45 


Črūss-immunē Individuāls 





0 5 10 15 20 25 30 
Time (dan) 


Fig 4. Shows the cross-immune individuals 
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shows the total time dependent population 


+ NSDF 
— — -Rka 
— — — ODE45 


Total population 





0 5 10 15 20 25 30 
Time {day} 


Fig 5. Shows the total time dependent population. 





Susceptible individuals 


0 5 10 15 20 25 30 
Time(day) 
Fig. 6. scr) versus t£ (solid line) MSGDTM, (dotted line) Runge-Kutta method. 
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250 


Infected individuals 
— p—s 
co „J 
S S 
—— 9 


Time(day) 
Fig. 7. 1(1) versus t (solid line) MSGDTM, (dotted line) Runge-Kutta method. 





Recovered individuals 


Time(day) 


Fig. 8. R(t) versus t (solid line) MSGDTM, (dotted line) Runge-Kutta method. 
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70 








Cross—immune individuals 
+) 
co 
OJ 


0 > 10 15 20 23 30 


Time(day) 


Fig. 9. ct) versus £ (solid line) MSGDTM, (dotted line) Runge-Kutta method. 


Total population 








Time(day) 


Fig. 10. y (+) versus £ (solid line) MSGDTM, (dotted line) Runge-Kutta method. 
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Susceptible individuals 





Time(day) 
Fig. 11. scr) versus t£ (solid line) a = 1.0 (dashed line) a = 0.95 , (dot-dashed line) a = 0.85. 


Infected individuals 





Time(day) 
Fig. 12. 7(t) versus £ (solid line) a =1.0 ( dashed line) a = 0.95 , (dot-dashed line) a = 0.85. 
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Recovered individuals 





Time(day) 
Fig.13. R(t) versus t£ (solid line) a =1.0 ( dashed line) a = 0.95 , (dot-dashed line) a = 0.85. 


Cross—immune individuals 





Time(day) 
Fig.14. ct) versus t (solid line) a = 1.0 ( dashed line) a = 0.95 , (dot-dashed line) a = 0.85. 
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Total population 





Time(day) 


Fig.15. w+) versus t (solid line) a =1.0 ( dashed line) a = 0.95 , (dot-dashed line) a = 0.85. 


6. Conclusion 


In this paper, a fractional order system for SIRC 
(Susceptible-Infected-Recovered-Cross-immune) epidemic model is studied and its 
approximate solution is presented using the multi-step generalized differential transform 
method (MSGDTM). 

The approximate solution obtained by multi-step generalized differential transform 
method are highly accurate and valid for a long time in the integer case. This method is 
very applicable and also this is a good approach for the solutions of differential equations 
of such order. 

This tool is the best one for modeling in science and engineering. 
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